# -*- coding: utf-8 -*-

"""
Comment:

"""
import numpy as np
from math import e
from scipy.stats import uniform, gumbel_r, gumbel_l
from gurobipy import *
from gurobipy import Model
from gurobipy import quicksum
import copy
from scipy.stats import binom
import time

def mnl(os2,prices,segQ,segments,segScale,segNP):
    prob_MNL = np.zeros([os2.shape[0],os2.shape[1]+1,segments.shape[0]])
    for l in segments:
        for i in np.arange(os2.shape[0]):
            prob_MNL[i,os2.shape[1],l] = 1/(e**(segNP[0,l])+np.sum(os2[i,:]*e**(1/segScale[l]*(segQ[:,l]-np.squeeze(prices[:])))))
            for j in np.arange(os2.shape[1]):
                prob_MNL[i,j,l] = os2[i,j]*(e**(1/segScale[l]*(segQ[j,l]-prices[j][0])))/(e**(segNP[0,l])+np.sum(os2[i,:]*e**(1/segScale[l]*(segQ[:,l]-np.squeeze(prices)))))
            prob_MNL[i,-1,l] = e**(segNP[0,l])/(e**(segNP[0,l])+np.sum(os2[i,:]*e**(1/segScale[l]*(segQ[:,l]-np.squeeze(prices)))))
    return prob_MNL

def getRandomstreams(products,segments,segLambdas,segQualityindices,segScale,timeperiods,runs,nEvalWSK,segNP,seed):
    nTimeperiods = timeperiods.size
    nProducts = products.size
    nSegments = segments.size
    nRuns = runs.size
    evals = np.arange(nEvalWSK)
    randomSegments = np.zeros((nTimeperiods,nRuns,nEvalWSK),dtype=np.int32)
    randomQualityindices = np.zeros((nTimeperiods,nProducts,nRuns,nEvalWSK))    
    randomNopurch = np.zeros((nTimeperiods,nRuns,nEvalWSK))
    
    U = uniform.rvs(loc=0,scale=1,size=[nTimeperiods,nRuns,nEvalWSK],random_state=seed) #Which segment?
    Q = np.zeros((nTimeperiods,nProducts+1,nRuns,nSegments,nEvalWSK)) #Which quality index
    for l in segments: #in for-loop because scale depends on segment
        randomNumbersOfSegment = gumbel_r.rvs(loc=0,scale=segScale[l],size=[nTimeperiods,nProducts+1,nRuns,nEvalWSK],random_state=seed)
        Q[:,:,:,l,:] = randomNumbersOfSegment
    for ev in evals:
        for run in runs:
            for t in timeperiods:
                segLambdasInPeriod = segLambdas #allgemeinfültig bei verschiedenen lambdas?!?!?!
                if U[t,run,ev] >= np.sum(segLambdasInPeriod):
                    randomSegments[t,run,ev] = -1000
                    randomQualityindices[t,:,run,ev] = -np.inf                
                else:
                    cumulativeProbs = np.cumsum(segLambdasInPeriod)
                    l = np.argmin(cumulativeProbs<U[t,run,ev])
                    randomSegments[t,run,ev] = l        
                    randomQualityindices[t,:,run,ev] = Q[t,:nProducts,run,l,ev] + segQualityindices[:,l]
                    randomNopurch[t,run,ev] = Q[t,nProducts,run,l,ev] + segNP[0,l]*segScale[l]
    return randomSegments,randomQualityindices,randomNopurch

def getRandomstreams_Buyer_j(nE,inst,runn,nTimeperiods,nRuns,nEvalWSK,c,segScale2,nProducts,j,offerset,l,prices,segQ,x):
    samplesize = 1
    nQSamples = samplesize*100000
    sampleQualityindices = np.zeros((samplesize,nProducts))   
    sampleNopurch = np.zeros((samplesize))
    t = 0
    while t < samplesize:
        Q = gumbel_r.rvs(loc=0,scale=segScale2,size=[nQSamples,nProducts+1],random_state=1+x)
        for tt in np.arange(nQSamples):
            q = Q[tt,:-1]
            utilities = q - np.squeeze(prices)                         #calculate new utilitites with replaced price
            utilities[offerset==0] = -np.inf                                  #set all products not in OS to -inf   
            y = np.argmax(utilities)                             
            if y == j:       
                sampleQualityindices[t,:] = Q[tt,:nProducts] + segQ[:,l]
                sampleNopurch[t] = Q[tt,nProducts]
                t += 1
            if t == samplesize:
                break
            else:
                next
        x +=1
    return sampleQualityindices,sampleNopurch,x

def cdlp(os,RS,capcon,cap,probMNL,nTimeperiods,products):
    ofs = np.arange(os.shape[0])
    mo = Model()
    ts = {}
    for o in np.arange(os.shape[0]):
        ts[o] = mo.addVar(lb=0,vtype='C',name="ts_"+str(o))
    obj = quicksum(np.squeeze(RS)[o]*ts[o] for o in ofs)
    mo.setObjective(obj)
    mo.update()
    #constraints
    constr_1 = {}
    constr_2 = {}
    constr_1 = mo.addConstr(quicksum(ts[o] for o in ofs) == nTimeperiods,name="constr_sumT")
    for h in np.arange(cap.shape[0]):
        constr_2[h] = mo.addConstr(quicksum(ts[o]*quicksum(capcon[p,h]*probMNL[o,p] for p in products) for o in ofs)<= cap[h,0],name="constr_kap_"+str(h))
    mo.update()
    #Further parameters and write()
    mo.setParam('OutputFlag',0)
    mo.setParam('Method',1)    
    mo.ModelSense = -1; #Maximization = -1
    mo.update()
    mo.optimize()
    mo.write('cdlp.lp')
    timeslots = np.zeros(os.shape[0])
    for o in ofs:
        timeslots[o] = ts[o].X
    oval = mo.objVal
    return timeslots, oval

def getRandomOSwithCDLP(nTimeperiods,nRuns,Ts,seed):
    mat = uniform.rvs(loc=0,scale=1,size=[nTimeperiods,nRuns],random_state=seed)
    OSNo = np.zeros(mat.shape)
    cumulativeProbs = np.cumsum(Ts/nTimeperiods)
    for run in np.arange(nRuns):
        for t in np.arange(nTimeperiods):
            l = np.argmin(cumulativeProbs<mat[t,run])
            OSNo[t,run] = l
    return OSNo

def efficientSets_largest_marginal_revenue_algorithm(Qt,RS,os):
    Qt = np.squeeze(Qt)
    RS = np.squeeze(RS)
    inSet = [0]
    effSets = list(inSet)
    posSets = np.squeeze(np.argwhere((Qt[:]>=Qt[inSet]) & (RS[:]>=RS[inSet])))
    posSets = posSets[posSets!=inSet]
    actSet = inSet
    while posSets.size != 0:
        actSet = np.argmax((RS[posSets]-RS[actSet])/(Qt[posSets]-Qt[actSet]))
        effSets.append(posSets[actSet])
        posSets = np.squeeze(np.argwhere((Qt[:]>=Qt[posSets[actSet]]) & (RS[:]>=RS[posSets[actSet]])))
        for i in effSets:
            posSets = posSets[posSets != i]
    return effSets

def calc_upsellpossibilities(uphierar,uniqSeg,counts2,nProducts,nRes,offerSets,prices,capcon,rem_cap):
    products = np.arange(nProducts)
    if nProducts == 2:
        nNewSegs = uniqSeg.shape[0]                                                    #number of new segments
        p_comp = np.zeros((nNewSegs,nProducts))
        i = 0
        uplist_M = np.zeros((nProducts,nProducts,nNewSegs))
        y_M = np.zeros((nProducts,nProducts,nNewSegs))
        uplist = []
        klist = []
        #segment specific prices - SINGLE-STEP CASCADING
        for nSeg in uniqSeg:
            u = []                                                 
            kli = []
            klis = 0
            #set price ranges
            for pro in products:
                if pro==0:       
                    p_comp[i,pro] = prices[nSeg[1].astype(int)]
                #prüfe auf Modulo, zudem höheres Produkt und erste Upsellmöglichkeit
                else:
                    p_comp[i,pro] = prices[pro][0]
                    u.append(pro)
                    kli.append(klis)
                    klis += 1
                    uplist_M[nSeg[1].astype(int),pro,i] = 1
                    y_M[nSeg[1].astype(int),pro,i] = counts2[i]
            i += 1
            uplist.append(u)
            klist.append(kli)    
    else:
        if uphierar == 'F2FnS2S':
            nNewSegs = uniqSeg.shape[0]                                                    #number of new segments
            p_comp = np.zeros((nNewSegs,nProducts))
            i = 0
            uplist_M = np.zeros((nProducts,nProducts,nNewSegs))
            y_M = np.zeros((nProducts,nProducts,nNewSegs))
            uplist = []
            klist = []
            #segment specific prices - SINGLE-STEP CASCADING
            for nSeg in uniqSeg:
                u = []                                                 
                kli = []
                klis = 0
                #set price ranges
                for pro in products:
                    if pro==nSeg[1].astype(int):       
                        p_comp[i,pro] = prices[nSeg[1].astype(int)]
                    #prüfe auf Modulo, zudem höheres Produkt und erste Upsellmöglichkeit
                    elif pro%2 == nSeg[1]%2 and pro>nSeg[1] and klis==0:
                        p_comp[i,pro] = prices[pro][0]
                        u.append(pro)
                        kli.append(klis)
                        klis += 1
                        uplist_M[nSeg[1].astype(int),pro,i] = 1
                        y_M[nSeg[1].astype(int),pro,i] = counts2[i]
                        next
                    #Ausschluss, wenn pro in höchster Klasse
                    elif (capcon[pro,:]*rem_cap[0])[nRes-1]==0 and np.where(capcon[pro,:]==1)[0]==nRes:
                        next
                    else:
                        p_comp[i,pro] = np.inf
                i += 1
                uplist.append(u)
                klist.append(kli)
            
        if uphierar == 'mixed':
            nNewSegs = uniqSeg.shape[0]                                                    #number of new segments
            p_comp = np.zeros((nNewSegs,nProducts))
            i = 0
            uplist_M = np.zeros((nProducts,nProducts,nNewSegs))
            y_M = np.zeros((nProducts,nProducts,nNewSegs))
            uplist = []
            klist = []
            #segment specific prices - SINGLE-STEP CASCADING (depending on offer set)
            for nSeg in uniqSeg:
                u = []                                                 
                kli = []
                klis = 0
                os = offerSets[nSeg[2]]
                #set price ranges
                for pro in products:
                    if pro==nSeg[1].astype(int):       
                        p_comp[i,pro] = prices[nSeg[1].astype(int)]
                    #beide Produkte in OS -> prüfe auf Modulo und in OS, zudem höheres Produkt und erste Upsellmöglichkeit
                    elif ((os@capcon*capcon[pro])[np.where(capcon[pro,:]==1)[0]]==2)[0] and pro%2 == nSeg[1].astype(int)%2 and os[pro]==1 and pro>nSeg[1].astype(int) and klis==0:
                        p_comp[i,pro] = prices[pro][0]
                        u.append(pro)
                        kli.append(klis)
                        klis += 1
                        uplist_M[nSeg[1].astype(int),pro,i] = 1
                        y_M[nSeg[1].astype(int),pro,i] = counts2[i]
                        next
                    #kein Produkt in OS -> prüfe auf Modulo, zudem höheres Produkt und erste Upsellmöglichkeit
                    elif ((os@capcon*capcon[pro])[np.where(capcon[pro,:]==1)[0]]==0)[0] and pro%2 == nSeg[1].astype(int)%2 and pro>nSeg[1].astype(int) and klis==0:
                        p_comp[i,pro] = prices[pro][0]
                        u.append(pro)
                        kli.append(klis)
                        klis += 1
                        uplist_M[nSeg[1].astype(int),pro,i] = 1
                        y_M[nSeg[1].astype(int),pro,i] = counts2[i]
                        next
                    #ein Produkt in OS -> Upsell dorthin
                    elif ((os@capcon*capcon[pro])[np.where(capcon[pro,:]==1)[0]]==1)[0] and os[pro]==1 and pro>nSeg[1].astype(int) and klis==0:
                        p_comp[i,pro] = prices[pro][0]
                        u.append(pro)
                        kli.append(klis)
                        klis += 1
                        uplist_M[nSeg[1].astype(int),pro,i] = 1
                        y_M[nSeg[1].astype(int),pro,i] = counts2[i]
                        next
                    #Ausschluss, wenn pro in höchster Klasse
                    elif (capcon[pro,:]*rem_cap[0])[nRes-1]==0 and np.where(capcon[pro,:]==1)[0]==nRes:
                        next
                    else:
                        p_comp[i,pro] = np.inf
                i += 1
                uplist.append(u)
                klist.append(kli)
    
    return uplist, klist, p_comp, uplist_M, y_M
            
def up_lin_mitPreisPs_mitSegmenten_Products(segScale,resources,uniqSeg,counts2,segQ,p_comp,cap,rem_cap,capcon,uplist,klist,offerSets,prices,nRes,method):
    
    model = Model("mip")
    
    #Param-Def
    nNewSegs = uniqSeg.shape[0]
    newSegs = np.arange(nNewSegs).astype(int)
    occ = np.squeeze(cap) - rem_cap                 #sold tickets
    p_comp = p_comp.astype(int) 
    
    plist = list()
    for l in newSegs:
        psub = []                                     #gurobi cannot handle min functions
        for k in uplist[l]:
            psub.append(list(np.linspace(p_comp[l,uniqSeg[l,1]],p_comp[l,k],(np.arange(p_comp[l,uniqSeg[l,1]],p_comp[l,k]).size)*1+1,endpoint=True)))        #1 steps
#            psub.append(list((prices[uniqSeg[l][1]][0],prices[uniqSeg[l][1]][0]+1)))                               #for small model to check constraints in lp.-file
        plist.append(psub)
    plist = [[[round(elem, 2) for elem in subsublist] for subsublist in sublist] for sublist in plist]
       
    comb_list = list()
    co_list = list()
    prob_list = list()
    
    for l in newSegs:
        comb =[]
        col = []
        prol = []
        os = offerSets[uniqSeg[l][2],:]
        sQ = segQ[:,uniqSeg[l][0]]
        sSc = segScale[uniqSeg[l][0]]
        for m in klist[l]:
            col.append(list(np.array(counts2[l]*(plist[l][m]-prices[uniqSeg[l,1]]))))
            prol.append(list(np.squeeze(np.array(1-((1+os@e**(1/sSc*(sQ-np.squeeze(prices))))/(1+os@e**(1/sSc*(sQ-np.squeeze(prices)))-os[np.where(prices==p_comp[l,uplist[l][m]])[0]]*e**(1/sSc*(sQ[np.where(prices==p_comp[l,uplist[l][m]])[0]]-p_comp[l,uplist[l][m]]))
                                        +e**(1/sSc*(sQ[np.where(prices==p_comp[l,uplist[l][m]])[0]]-plist[l][m]))))))))
            comb.append(list(np.squeeze(np.array(counts2[l]*(plist[l][m]-prices[uniqSeg[l,1]]))*np.array(1-((1+os@e**(1/sSc*(sQ-np.squeeze(prices))))/(1+os@e**(1/sSc*(sQ-np.squeeze(prices)))-os[np.where(prices==p_comp[l,uplist[l][m]])[0]]*e**(1/sSc*(sQ[np.where(prices==p_comp[l,uplist[l][m]])[0]]-p_comp[l,uplist[l][m]]))
                                        +e**(1/sSc*(sQ[np.where(prices==p_comp[l,uplist[l][m]])[0]]-plist[l][m]))))))))
        comb_list.append(comb)
        co_list.append(col)
        prob_list.append(prol)
        
    #Var-Def and Tuple for OF
    xkpl = {}
    tup = list()
    for l in newSegs:
        for m in klist[l]:
            for p in plist[l][m]:
                xkpl[m,p,l] = model.addVar(lb=0,ub=1,vtype='B',name="xkpl_"+str(m)+str(p)+str(l))
                tup.append((m,p,l))
    model.update() 
    
    #Objective Function
    obj = quicksum(comb_list[l][m][int(p-prices[uniqSeg[l][1]])*1]*xkpl[m,p,l] for (m,p,l) in tup)
    model.setObjective(obj)
    model.update() 
    
    #Constraints
    constr_1 = {}
    constr_2 = {}
    for l in newSegs:              #maximal ein Upsellangebot pro Segment
        constr_1[l] = model.addConstr(quicksum(quicksum(xkpl[m,p,l] for p in plist[l][m]) for m in klist[l])<=1,name="constr1_"+str(l))
    
    counts3 = np.zeros((uniqSeg.shape[0],resources.shape[0]))
    counts4 = copy.deepcopy(counts3)
    n = np.zeros(uniqSeg.shape[0])
    newSegs_res = np.zeros((uniqSeg.shape[0],resources.shape[0]))
    for l in newSegs:
        counts3[l,np.array(klist[l])+np.where(capcon[uniqSeg[l,1],:]==1)] = counts2[l]
        counts4[l,np.where(capcon[uniqSeg[l,1],:]==1)] = counts2[l]
        n[l] = np.where(capcon[uniqSeg[l,1],:]==1)[0]
        newSegs_res[l,np.array(klist[l])+np.where(capcon[uniqSeg[l,1],:]==1)[0]] = newSegs[l]
    nSR = list()
    n = n.astype(int)
    for h in resources:
        nSR.append(list(newSegs_res[:,h]))
    newSegsRes = [[int(elem) for elem in sublist if elem != 0] for sublist in nSR]
    if rem_cap[0][nRes-1]==0:
            h_scope = np.arange(nRes)[1:nRes-1]
    else:
        h_scope = np.arange(nRes)[1:]
    
    if method == 1:
        [newSegsRes[np.squeeze(np.where(capcon[uniqSeg[0,1]]==1))].insert(0,0)]
        for h in h_scope:
            constr_2[h] = model.addConstr(occ[0][h]                                                                                           #momentane Auslastung
                    + quicksum(quicksum(counts3[l,d]*quicksum(prob_list[l][d-n[l]][int((p-prices[uniqSeg[l][1]])*1)]*xkpl[d-n[l],p,l] for p in plist[l][d-n[l]])for l in newSegsRes[d])for d in list([h-1]))    #Zugaenge
                    - quicksum(quicksum(counts4[l,h]*quicksum(prob_list[l][m][int((p-prices[uniqSeg[l][1]])*1)]*xkpl[m,p,l] for p in plist[l][m])for m in klist[l])for l in newSegs)      #Abgaenge
                    <= cap[h][0],name="constr2_"+str(h))
    elif method == 2:
        [newSegsRes[h].insert(0,0) for h in resources[:-1]]
        for h in h_scope:
            constr_2[h] = model.addConstr(occ[0][h]                                                                                           #momentane Auslastung
                    + quicksum(quicksum(counts3[l,d]*quicksum(prob_list[l][d-n[l]][int((p-prices[uniqSeg[l][1]])*1)]*xkpl[d-n[l],p,l] for p in plist[l][d-n[l]])for l in newSegsRes[d])for d in list(range(h)))    #Zugaenge
                    - quicksum(quicksum(counts4[l,h]*quicksum(prob_list[l][m][int((p-prices[uniqSeg[l][1]])*1)]*xkpl[m,p,l] for p in plist[l][m])for m in klist[l])for l in newSegs)      #Abgaenge
                    <= cap[h][0],name="constr2_"+str(h))
                
    model.update()
    model.setParam('OutputFlag',0)
    model.setParam('Method',1)
    model.setParam('MIPGap',0.01) 
    model.ModelSense = -1; #Maximization = -1
#    model.setParam('NonConvex',2)
    model.update()
    time_start = time.perf_counter_ns()
    model.optimize()
    tsub = (time.perf_counter_ns() - time_start)
    count = newSegs.shape[0]
    values = np.zeros((count,4))
    i=0
    for l in newSegs:
        for m in klist[l]:
            for p in plist[l][m]:
                if xkpl[m,p,l].x == 1:
                    values[i,:] = np.array([l,uplist[l][m],p,xkpl[m,p,l].x])
                    i += 1
    values = values[0:i]
    revenue = model.objVal
    model.write("upprob.lp")
    return values, revenue, tsub

def upsellproblem_relaxed_CapConst(segScale,resources,uniqSeg,counts2,segQ,p_comp,cap,rem_cap,capcon,uplist,klist,offerSets,prices,nRes,weight,method):
    
    model = Model("mip")
    
    #Param-Def
    nNewSegs = uniqSeg.shape[0]
    newSegs = np.arange(nNewSegs).astype(int)
    occ = np.squeeze(cap) - rem_cap                 #sold tickets
    p_comp = p_comp.astype(int) 
    
    plist = list()
    for l in newSegs:
        psub = []                                     #gurobi cannot handle min functions
        for k in uplist[l]:
            psub.append(list(np.linspace(p_comp[l,uniqSeg[l,1]],p_comp[l,k],(np.arange(p_comp[l,uniqSeg[l,1]],p_comp[l,k]).size)*1+1,endpoint=True)))        #1 steps
#            psub.append(list((prices[uniqSeg[l][1]][0],prices[uniqSeg[l][1]][0]+1)))                               #for small model to check constraints in lp.-file
        plist.append(psub)
    plist = [[[round(elem, 2) for elem in subsublist] for subsublist in sublist] for sublist in plist]
       
    comb_list = list()
    co_list = list()
    prob_list = list()
    
    for l in newSegs:
        comb =[]
        col = []
        prol = []
        os = offerSets[uniqSeg[l][2],:]
        sQ = segQ[:,uniqSeg[l][0]]
        sSc = segScale[uniqSeg[l][0]]
        for m in klist[l]:
            col.append(list(np.array(counts2[l]*(plist[l][m]-prices[uniqSeg[l,1]]))))
            prol.append(list(np.squeeze(np.array(1-((1+os@e**(1/sSc*(sQ-np.squeeze(prices))))/(1+os@e**(1/sSc*(sQ-np.squeeze(prices)))-os[np.where(prices==p_comp[l,uplist[l][m]])[0]]*e**(1/sSc*(sQ[np.where(prices==p_comp[l,uplist[l][m]])[0]]-p_comp[l,uplist[l][m]]))
                                        +e**(1/sSc*(sQ[np.where(prices==p_comp[l,uplist[l][m]])[0]]-plist[l][m]))))))))
            comb.append(list(np.squeeze(np.array(counts2[l]*(plist[l][m]-prices[uniqSeg[l,1]]))*np.array(1-((1+os@e**(1/sSc*(sQ-np.squeeze(prices))))/(1+os@e**(1/sSc*(sQ-np.squeeze(prices)))-os[np.where(prices==p_comp[l,uplist[l][m]])[0]]*e**(1/sSc*(sQ[np.where(prices==p_comp[l,uplist[l][m]])[0]]-p_comp[l,uplist[l][m]]))
                                        +e**(1/sSc*(sQ[np.where(prices==p_comp[l,uplist[l][m]])[0]]-plist[l][m]))))))))
        comb_list.append(comb)
        co_list.append(col)
        prob_list.append(prol)
    
    #Variables/Calculation of sets for Constraints
    counts3 = np.zeros((uniqSeg.shape[0],resources.shape[0]))
    counts4 = copy.deepcopy(counts3)
    n = np.zeros(uniqSeg.shape[0])
    newSegs_res = np.zeros((uniqSeg.shape[0],resources.shape[0]))
    for l in newSegs:
        counts3[l,np.array(klist[l])+np.where(capcon[uniqSeg[l,1],:]==1)] = counts2[l]
        counts4[l,np.where(capcon[uniqSeg[l,1],:]==1)] = counts2[l]
        n[l] = np.where(capcon[uniqSeg[l,1],:]==1)[0]
        newSegs_res[l,np.array(klist[l])+np.where(capcon[uniqSeg[l,1],:]==1)[0]] = newSegs[l]
    nSR = list()
    n = n.astype(int)
    for h in resources:
        nSR.append(list(newSegs_res[:,h]))
    newSegsRes = [[int(elem) for elem in sublist if elem != 0] for sublist in nSR]
    if rem_cap[0][nRes-1]==0:
            h_scope = np.arange(nRes)[1:nRes-1]
    else:
        h_scope = np.arange(nRes)[1:]
        
    #Var-Def and Tuple for OF
    xkpl = {}
    punish = {}
    tup = list()
    for l in newSegs:
        for m in klist[l]:
            for p in plist[l][m]:
                xkpl[m,p,l] = model.addVar(lb=0,ub=1,vtype='B',name="xkpl_"+str(m)+str(p)+str(l))
                tup.append((m,p,l))
    for h in h_scope:
        punish[h] = model.addVar(lb=0,vtype='C',name="punish_"+str(h))
    model.update() 
    
    #Objective Function
    obj = quicksum(comb_list[l][m][int(p-prices[uniqSeg[l][1]])*1]*xkpl[m,p,l] for (m,p,l) in tup) - quicksum(weight*punish[h] for h in h_scope)
    model.setObjective(obj)
    model.update() 
    
    #Constraints
    constr_1 = {}
    constr_2 = {}
    for l in newSegs:              #maximal ein Upsellangebot pro Segment
        constr_1[l] = model.addConstr(quicksum(quicksum(xkpl[m,p,l] for p in plist[l][m]) for m in klist[l])<=1,name="constr1_"+str(l))
    
    if method == 1:
        [newSegsRes[np.squeeze(np.where(capcon[uniqSeg[0,1]]==1))].insert(0,0)]
        for h in h_scope:
            constr_2[h] = model.addConstr(occ[0][h]                                                                                           #momentane Auslastung
                    + quicksum(quicksum(counts3[l,d]*quicksum(prob_list[l][d-n[l]][int((p-prices[uniqSeg[l][1]])*1)]*xkpl[d-n[l],p,l] for p in plist[l][d-n[l]])for l in newSegsRes[d])for d in list([h-1]))    #Zugaenge
                    - quicksum(quicksum(counts4[l,h]*quicksum(prob_list[l][m][int((p-prices[uniqSeg[l][1]])*1)]*xkpl[m,p,l] for p in plist[l][m])for m in klist[l])for l in newSegs)      #Abgaenge
                    <= cap[h][0]+punish[h],name="constr2_"+str(h))
    elif method == 2:
        [newSegsRes[h].insert(0,0) for h in resources[:-1]]
        for h in h_scope:
            constr_2[h] = model.addConstr(occ[0][h]                                                                                           #momentane Auslastung
                    + quicksum(quicksum(counts3[l,d]*quicksum(prob_list[l][d-n[l]][int((p-prices[uniqSeg[l][1]])*1)]*xkpl[d-n[l],p,l] for p in plist[l][d-n[l]])for l in newSegsRes[d])for d in list(range(h)))    #Zugaenge
                    - quicksum(quicksum(counts4[l,h]*quicksum(prob_list[l][m][int((p-prices[uniqSeg[l][1]])*1)]*xkpl[m,p,l] for p in plist[l][m])for m in klist[l])for l in newSegs)      #Abgaenge
                    <= cap[h][0],name="constr2_"+str(h))
                
    model.update()
    model.setParam('OutputFlag',0)
    model.setParam('Method',1)
    model.setParam('MIPGap',0.01) 
    model.ModelSense = -1; #Maximization = -1
#    model.setParam('NonConvex',2)
    model.update()
    time_start = time.perf_counter_ns()
    model.optimize()
    tsub = (time.perf_counter_ns() - time_start)
    count = newSegs.shape[0]
    values = np.zeros((count,4))
    pp = np.zeros(resources.size)
    i=0
    for l in newSegs:
        for m in klist[l]:
            for p in plist[l][m]:
                if xkpl[m,p,l].x == 1:
                    values[i,:] = np.array([l,uplist[l][m],p,xkpl[m,p,l].x])
                    i += 1
    for h in h_scope:
        pp[h] = punish[h].x
    values = values[0:i]
    revenue = model.objVal
    model.write("upprob.lp")
    return values, revenue, pp, tsub

def up_lin_mean(segScale,resources,uniqSeg,counts2,segQ,p_comp,cap,rem_cap,capcon,uplist,klist,offerSets,prices,nRes,method):
    
    model = Model("mip")
    
    #Param-Def
    nNewSegs = uniqSeg.shape[0]
    newSegs = np.arange(nNewSegs).astype(int)
    occ = np.squeeze(cap) - rem_cap                 #sold tickets
    p_comp = p_comp.astype(int) 
    
    plist = list()
    for l in newSegs:
        psub = []                                     #gurobi cannot handle min functions
        for k in uplist[l]:
            psub.append(list(np.linspace(p_comp[l,uniqSeg[l,1]],p_comp[l,k],(np.arange(p_comp[l,uniqSeg[l,1]],p_comp[l,k]).size)*1+1,endpoint=True)))        #1 steps
#            psub.append(list((prices[uniqSeg[l][1]][0],prices[uniqSeg[l][1]][0]+1)))                               #for small model to check constraints in lp.-file
        plist.append(psub)
    plist = [[[round(elem, 2) for elem in subsublist] for subsublist in sublist] for sublist in plist]
       
    comb_list = list()
    co_list = list()
    prob_list = list()
    
    for l in newSegs:
        comb =[]
        col = []
        prol = []
        os = offerSets[uniqSeg[l][2],:]
        sQ = segQ[:,uniqSeg[l][0]]
        sSc = segScale[uniqSeg[l][0]]
        for m in klist[l]:
            col.append(list(np.array(counts2[l]*(plist[l][m]-prices[uniqSeg[l,1]]))))
            prol.append(list(np.squeeze(np.array(1-((1+os@e**(1/sSc*(sQ-np.squeeze(prices))))/(1+os@e**(1/sSc*(sQ-np.squeeze(prices)))-os[np.where(prices==p_comp[l,uplist[l][m]])[0]]*e**(1/sSc*(sQ[np.where(prices==p_comp[l,uplist[l][m]])[0]]-p_comp[l,uplist[l][m]]))
                                        +e**(1/sSc*(sQ[np.where(prices==p_comp[l,uplist[l][m]])[0]]-plist[l][m]))))))))
            comb.append(list(np.squeeze(np.array(counts2[l]*(plist[l][m]-prices[uniqSeg[l,1]]))*np.array(1-((1+os@e**(1/sSc*(sQ-np.squeeze(prices))))/(1+os@e**(1/sSc*(sQ-np.squeeze(prices)))-os[np.where(prices==p_comp[l,uplist[l][m]])[0]]*e**(1/sSc*(sQ[np.where(prices==p_comp[l,uplist[l][m]])[0]]-p_comp[l,uplist[l][m]]))
                                        +e**(1/sSc*(sQ[np.where(prices==p_comp[l,uplist[l][m]])[0]]-plist[l][m]))))))))
        comb_list.append(comb)
        co_list.append(col)
        prob_list.append(prol)
        
    #Var-Def and Tuple for OF
    xkpl = {}
    tup = list()
    for l in newSegs:
        for m in klist[l]:
            for p in plist[l][m]:
                xkpl[m,p,l] = model.addVar(lb=0,ub=1,vtype='B',name="xkpl_"+str(m)+str(p)+str(l))
                tup.append((m,p,l))
    model.update() 
    
    #Objective Function
    obj = quicksum(comb_list[l][m][int(p-prices[uniqSeg[l][1]])*1]*xkpl[m,p,l] for (m,p,l) in tup)
    model.setObjective(obj)
    model.update() 
    
    #Constraints
    constr_1 = {}
    constr_2 = {}
    for l in newSegs:              #maximal ein Upsellangebot pro Segment
        constr_1[l] = model.addConstr(quicksum(quicksum(xkpl[m,p,l] for p in plist[l][m]) for m in klist[l])<=1,name="constr1_"+str(l))
    
    counts3 = np.zeros((uniqSeg.shape[0],resources.shape[0]))
    counts4 = copy.deepcopy(counts3)
    n = np.zeros(uniqSeg.shape[0])
    newSegs_res = np.zeros((uniqSeg.shape[0],resources.shape[0]))
    for l in newSegs:
        counts3[l,np.array(klist[l])+np.where(capcon[uniqSeg[l,1],:]==1)] = counts2[l]
        counts4[l,np.where(capcon[uniqSeg[l,1],:]==1)] = counts2[l]
        n[l] = np.where(capcon[uniqSeg[l,1],:]==1)[0]
        newSegs_res[l,np.array(klist[l])+np.where(capcon[uniqSeg[l,1],:]==1)[0]] = newSegs[l]
    nSR = list()
    n = n.astype(int)
    for h in resources:
        nSR.append(list(newSegs_res[:,h]))
    newSegsRes = [[int(elem) for elem in sublist if elem != 0] for sublist in nSR]
    if rem_cap[0][nRes-1]==0:
            h_scope = np.arange(nRes)[1:nRes-1]
    else:
        h_scope = np.arange(nRes)[1:]
    
    if method == 1:
        [newSegsRes[np.squeeze(np.where(capcon[uniqSeg[0,1]]==1))].insert(0,0)]
        for h in h_scope:
            constr_2[h] = model.addConstr(occ[0][h]                                                                                           #momentane Auslastung
                    + quicksum(quicksum(counts3[l,d]*quicksum(prob_list[l][d-n[l]][int((p-prices[uniqSeg[l][1]])*1)]*xkpl[d-n[l],p,l] for p in plist[l][d-n[l]])for l in newSegsRes[d])for d in list([h-1]))    #Zugaenge
                    - quicksum(quicksum(counts4[l,h]*quicksum(prob_list[l][m][int((p-prices[uniqSeg[l][1]])*1)]*xkpl[m,p,l] for p in plist[l][m])for m in klist[l])for l in newSegs)      #Abgaenge
                    <= cap[h][0],name="constr2_"+str(h))
    elif method == 2:
        [newSegsRes[h].insert(0,0) for h in resources[:-1]]
        for h in h_scope:
            constr_2[h] = model.addConstr(occ[0][h]                                                                                           #momentane Auslastung
                    + quicksum(quicksum(counts3[l,d]*quicksum(prob_list[l][d-n[l]][int((p-prices[uniqSeg[l][1]])*1)]*xkpl[d-n[l],p,l] for p in plist[l][d-n[l]])for l in newSegsRes[d])for d in list(range(h)))    #Zugaenge
                    - quicksum(quicksum(counts4[l,h]*quicksum(prob_list[l][m][int((p-prices[uniqSeg[l][1]])*1)]*xkpl[m,p,l] for p in plist[l][m])for m in klist[l])for l in newSegs)      #Abgaenge
                    <= cap[h][0],name="constr2_"+str(h))
                
    model.update()
    model.setParam('OutputFlag',0)
    model.setParam('Method',1)
    model.setParam('MIPGap',0.01) 
    model.ModelSense = -1; #Maximization = -1
#    model.setParam('NonConvex',2)
    model.update()
    time_start = time.perf_counter_ns()
    model.optimize()
    tsub = (time.perf_counter_ns() - time_start)
    count = newSegs.shape[0]
    values = np.zeros((count,4))
    i=0
    for l in newSegs:
        for m in klist[l]:
            for p in plist[l][m]:
                if xkpl[m,p,l].x == 1:
                    values[i,:] = np.array([l,uplist[l][m],p,xkpl[m,p,l].x])
                    i += 1
    values = values[0:i]
    revenue = model.objVal
    model.write("upprob.lp")
    return values, revenue, tsub